---
title: "Vignette"
author: "Dai Yamao"
date: "`r format(Sys.time(), '%Y-%m-%d')`"
output:
  html_document: default
  pdf_document: default
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE, fig.width = 10, fig.height = 5)
rm(list=ls())

require(memisc)
require(stringi)
library(stargazer)
library(coefplot)
library(interplot)
library(sjPlot)
library(sjmisc)
library(sjlabelled)
library(stargazer)
library(stringr)
library(tidyverse)
library(ggplot2)
library(psych)
library(margins)
library(dplyr)
library(modelsummary)
library(makedummies)
library(ggeffects)
library(estimatr)
library(ggpubr)

dat <- read_csv("data/vignette_data.csv")

dat <- dat %>% 
  mutate_at(c("Riot_Initical", "Riot_Political", "Riot_Iran", 
              "Remain_Initical", "Remain_Political", "Remain_Iran"), ~ replace(., is.na(.), 0))

dat <- dat %>% 
  mutate(sunni = ifelse(Sect == 1, 1, 0),
         shia = ifelse(Sect == 2, 1, 0),
         kurd = ifelse(Sect == 4, 1, 0),
         christian = ifelse(Sect == 3, 1, 0),
         minority = ifelse(Sect == 5, 1, 0))

dat$Q3_1[dat$Q3_1 == 6] <- NA
dat$Q3_2[dat$Q3_2 == 6] <- NA
dat$Q3_3[dat$Q3_3 == 6] <- NA
dat$Q3_7[dat$Q3_7 == 6] <- NA

dat <- dat %>% 
  mutate(trust_president = 6 - Q3_1,
         trust_pm = 6 - Q3_2,
         trust_parliament = 6 - Q3_3,
         trust_party = 6 - Q3_7)
dat <- dat %>% 
  mutate(trust_gov = (trust_president + trust_pm + trust_parliament + trust_party) /4)

dat <- dat %>% 
  mutate(Peace_protest = case_when(Remain_Initical == 1 | Remain_Political == 1 | Remain_Iran == 1 ~ 1,
                                   Riot_Initical == 1 | Riot_Political == 1 | Riot_Iran == 1 ~ 0))
dat$Peace_protest <- as.factor(dat$Peace_protest)

dat <- dat %>% 
  mutate(Internal_political_reform = case_when(Remain_Initical == 1 | Remain_Political == 1 | Riot_Initical == 1 |  Riot_Political == 1 ~ 1,
                                               Remain_Iran == 1 | Riot_Iran == 1 ~ 0))
dat$Internal_political_reform <- as.factor(dat$Internal_political_reform)

dat <- dat |>
  mutate(Treatment = case_when(Remain_Initical == 1 ~ "Peace_initial",
                               Remain_Political == 1 ~ "Peace_internal",
                               Remain_Iran == 1 ~ "Peace_external",
                               Riot_Initical == 1 ~ "Riot_initial",
                               Riot_Political == 1 ~ "Riot_internal",
                               Riot_Iran == 1 ~ "Riot_external"))
dat$Treatment <- as.factor(dat$Treatment) 
```

# overview
```{r}
stargazer(as.data.frame(dat), type = "text")
dat[,c(41:46)] <- data.frame(lapply(dat[,c(41:46)],as.factor))
str(dat)

dat %>% 
  filter(!is.na(sympathy)) %>% 
  ggplot() +
  geom_histogram(aes(x = sympathy), color = "white", bins = 4) +
  labs(x = "Sympathy toward protests", y = "Frequency")
```

# regression
```{r}
reg01 <- lm(formula = sympathy ~ Riot_Initical + Riot_Political +  
              Remain_Initical + Remain_Political + Remain_Iran,
            data = dat)
summary(reg01)

reg02 <- lm(formula = sympathy ~ Riot_Initical + Riot_Political +  
              Remain_Initical + Remain_Political + Remain_Iran +
              gender + age + Income + Education,
            data = dat)
summary(reg02)
stargazer(reg02, type = "text")

reg03 <- lm(formula = sympathy ~ Riot_Initical + Riot_Political +  
              Remain_Initical + Remain_Political + Remain_Iran +
              shia + sunni + kurd + christian +
              gender + age + Income + Education,
            data = dat)
summary(reg03)
stargazer(reg03, type = "text")
```

# regression: subset
```{r}
# sunni
reg04 <- lm(formula = sympathy ~ Riot_Initical + Riot_Political +  
              Remain_Initical + Remain_Political + Remain_Iran +
              gender + age + Income + Education,
            subset(dat, Sect == 1))
summary(reg04)

# shia
reg05 <- lm(formula = sympathy ~ Riot_Initical + Riot_Political +  
              Remain_Initical + Remain_Political + Remain_Iran +
              gender + age + Income + Education,
            subset(dat, Sect == 2))
summary(reg05)

# kurd
reg06 <- lm(formula = sympathy ~ Riot_Initical + Riot_Political +  
              Remain_Initical + Remain_Political + Remain_Iran +
              gender + age + Income + Education,
            subset(dat, Sect == 4))
summary(reg06)

```

# peace and internal aggregarion
```{r}
reg07 <- lm(formula = sympathy ~ Peace_protest +
              shia + sunni + kurd + christian +
              gender + age + Income + Education,
            data = dat)
summary(reg07)
stargazer(reg07, type = "text")

reg08 <- lm(formula = sympathy ~ Internal_political_reform +
              shia + sunni + kurd + christian +
              gender + age + Income + Education,
            data = dat)
summary(reg08)
stargazer(reg08, type = "text")

reg20 <- lm(formula = sympathy ~ Treatment +
              shia + sunni + kurd + christian +
              gender + age + Income + Education,
            data = dat)
summary(reg20)
```

# regression: table
```{r}
reg_table <- mtable("Model 1" = reg03, "Model 2(Peace Protest)" = reg07, 
                    "Model 3(Internal issue)" = reg08,
               summary.stats=c("adj. R-squared","F","p", "N"))
print(reg_table)
write_html(reg_table, file = "result/result_vignette.html")
```

# regression: trust
```{r}
reg09 <- lm(formula = sympathy ~ Riot_Initical + Riot_Political +  
              Remain_Initical + Remain_Political + Remain_Iran +
              trust_president + trust_pm + trust_parliament + trust_party +
              gender + age + Income + Education,
            data = dat)
summary(reg09)

reg10 <- lm(formula = sympathy ~ Riot_Initical + Riot_Political +  
              Remain_Initical + Remain_Political + Remain_Iran +
              trust_gov +
              gender + age + Income + Education,
            data = dat)
summary(reg10)
stargazer(reg10, type = "text")
```

# marginal effect
```{r}
plot_model(reg01, type = "pred", terms = c("Remain_Political"),
           axis.labels = "Treatment", title = "")
plot_model(reg03, type = "pred", terms = c("Remain_Political"),
           axis.labels = "Treatment", title = "")


plot_model(reg03, type = "pred", terms = c("Remain_Political", "gender", "Education"),
           axis.labels = "Treatment", title = "")
plot_model(reg03, type = "pred", terms = c("Remain_Political", "gender","Income"),
           axis.labels = "Treatment", title = "")
plot_model(reg03, type = "pred", terms = c("Education", "Remain_Political"),
           axis.labels = "Treatment", title = "")
plot_model(reg03, type = "pred", terms = c("Income", "Remain_Political"),
           axis.labels = "Treatment", title = "")

plot_model(reg09, type = "pred", terms = c("trust_party", "Remain_Political"),
           axis.labels = "Treatment", title = "")
plot_model(reg09, type = "pred", terms = c("Remain_Political", "trust_party"),
           axis.labels = "Treatment", title = "")
plot_model(reg10, type = "pred", terms = c("trust_gov", "Remain_Political"),
           axis.labels = "Treatment", title = "")

plot_model(reg20, type = "pred", terms = c("Treatment"),
           axis.labels = "Treatment", title = "")
```

# marginal effect
```{r}
p1 <- plot_model(reg03, type = "pred", terms = c("Remain_Political"),
           axis.labels = "Treatment", title = "")
p2 <- plot_model(reg03, type = "pred", terms = c("Remain_Initical"),
           axis.labels = "Treatment", title = "")
p3 <- plot_model(reg03, type = "pred", terms = c("Remain_Iran"),
           axis.labels = "Treatment", title = "")
p4 <- plot_model(reg03, type = "pred", terms = c("Riot_Political"),
           axis.labels = "Treatment", title = "")
p5 <- plot_model(reg03, type = "pred", terms = c("Riot_Initical"),
           axis.labels = "Treatment", title = "")

ggarrange(p5, 
          p4 + theme(axis.title.y = element_blank()), 
          p2 + theme(axis.title.y = element_blank()), 
          p1, 
          p3 + theme(axis.title.y = element_blank()),
          nrow = 2, ncol = 3, align = "hv")
```

# marginal effect
```{r}
p6 <- plot_model(reg07, type = "pred", terms = c("Peace_protest"),
           axis.labels = "Treatment", title = "")
p7 <- plot_model(reg08, type = "pred", terms = c("Internal_political_reform"),
           axis.labels = "Treatment", title = "")

ggarrange(p6, 
          p7 + theme(axis.title.y = element_blank()), 
          nrow = 1, ncol = 2, align = "hv")
```
